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Abstract 

In this paper, we provide the 0(e) corrections to the hydrodynamic model de- 
rived by Degond & Motsch from a kinetic version of the model by Vicsek &: coauthors 
describing flocking biological agents. The parameter e stands for the ratio of the 
microscopic to the macroscopic scales. The 0{e) corrected model involves diffusion 
terms in both the mass and velocity equations as well as terms which are quadratic 
functions of the first order derivatives of the density and velocity. The derivation 
method is based on the standard Chapman-Enskog theory, but is significantly more 
complex than usual due to both the non-isotropy of the fluid and the lack of mo- 
mentum conservation. 
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1 Introduction 



This paper is a development of a previous work [261 123 about continuum models of self- 
propelled particles subject to alignment interaction. This class of models describes swarm- 
ing behaviour among biological species and attempts at providing a simplified theoretical 
framework to experimental observations (see recent observations in Refs. |3l El IH1 l32l HO] ). 

The starting point is a particle model (or Individual-Based Model (IBM)) discussed 
in e.g. Refs. [2l |T6l |39l HI], where the interaction between biological agents such as fish 
or birds is described by three interaction ranges: a close range where repulsion occurs, 
a long range where attraction prevails and a medium range where agents tend to align 
with each other. In the Vicsek model, [53] the alignment interaction is singled out and 
analyzed. More precisely, each particle moves with a constant and uniform speed and 
aligns with the average direction of all neighbours within an interaction distance i?, up 
to some angular fluctuation. Vicsek and co-authors [53] show that phase transitions from 
disorder to order appear as the noise intensity decreases or the density increases. This 
model has triggered a wealth of publications [H HI [131 ESI ESI [3l] and given rise to various 
variants [m 150] . 

The task of deriving kinetic (Boltzmann-like) or continuum (fluid-like) models from 
this model has been undertaken by various approaches (see Refs. [6l |16] for kinetic models 
and Refs. [6l [2T1 |23lll6lll7] for fluid models). However, these models are based on physical 
arguments and the mathematical approach of Ref. [26l [27| leads to a different type of 
model: indeed, this model is not of diffusive nature like in Refs. [HI EH |23l HH] and is 
local, by contrast to Ref. [17]. The differences arise because the model of Refs. [261 EZ] is 
derived in the large time and space scale limit in which, at leading order, the interactions 
are local and the diffusivities, negligible. 

The goal of this paper is to investigate how the model of Ref. [26l [27] must be 
adapted to account for small but finite nonlocality and diffusivity. We will see that the 
resulting model involves complex effects due to the non-isotropy of the fluid (by contrast 
to standard fluids). This point will be further developed in section |2l Similar to the 
rarefied gas dynamics case (see e.g. Refs. [9l [121 1211 US]), the derivation is based on 
the Chapman-Enskog expansion method. However, in the present case, the computations 
are significantly more complex because of the non-isotropy of the fluid and of the lack of 
momentum conservation. 

The model derived in Ref. [261 l27] has been extended in Ref. [31] to account for 
anisotropic vision and density- dependent interaction frequency. Its numerical resolution 
is performed in Ref. [44j. 

Other Individual-Based Models involving attraction-repulsion interactions can be found 
in Refs. [2S1 ISOliSlilS] and continuum models, in Refs. [71[ini[I51l2nil251lll[5IllS2]. 
The existence of flocking or non-flocking behaviour for the Cucker-Smale model[T7 1 [T8 l [19] 
(which is similar to the Vicsek model but without noise nor speed constraint) has received 
a great deal of attention[ni |29l [351 |36l [371 |38lll8] . 



2 



2 Position of the problem and main result 



In Refs. [26l |22] the following continuum model of self-driven particles with alignment 
interaction has been derived: 



+ V ■ (cipfi) = 0, (2.1) 
p {dtn + C2(fi ■ V)fi) + C3 (Id - 1] ® n)Vp = 0, (2.2) 



where p = p{x, t) is the number density of the particles and Q{x, t) is the direction of their 
average velocity, which satisfies \Vt\ = 1. ci, C2 and C3 are constants which are computed 
from the underlying microscopic dynamics, where C2 < Ci < 1 and Ci, C2, C3 > 0. The 
matrix (Id — i7 (g) i7) denotes the orthogonal projection matrix onto the plane orthogonal 
to Q. The notations Id and ® respectively refer to the identity matrix and to the tensor 
product. 

The derivation of this model, according to Ref. [27], proceeds through several (formal) 
asymptotic limits. The starting point is a time-discrete particle model proposed by Vicsek 
and CO- authors [53]. In the Vicsek model, the particles move with a constant speed and, 
at discrete times, align their velocities to the mean velocity of their neighbours, up to 
some small noise. In Ref. [27], a continuous in time version of this particle model is first 
proposed in which the alignment interaction is modeled through a relaxation term and 
the noise, by a Brownian motion on the particle velocities. The formal mean-field limit 
of this time-continuous particle model (as the number of particles tends to infinity) leads 
to the following nonlinear Fokker-Planck model: 

dJ + u-Vf + V^- {Ff) = dA^f, (2.3) 
F(x, w, t) = z/(cos 6) {Id — u ® uj)tu{x, t), (2.4) 

ui{x,u,t) = I'^i^'^L J{x,t)=f K{\x -y\)v f{y,v,t)dydv . (2.5) 

Here f{x, 00, t) is the particle distribution function depending on the space variable x G M'^, 
the velocity direction a; G §^ and the time t. dis a scaled diffusion constant associated with 
the Brownian noise, and F{x, u, t) is the mean-field interaction force between the particles 
which depends on an interaction frequency u. This force tends to align the particles to the 
direction u which is the direction of the particle flux J in a neighbourhood of x weighted 
by the kernel K. Typically, if K is the indicator function of the ball centered at and 
of radius R then J is the particle flux integrated over a ball centered at x and of radius 
R. The matrix (Id — u ® u) is the projection matrix onto the normal plane to u and we 
assume that the collision frequency may depend on cos^ = {u ■ u), i.e. on the cosine of 
the angle between u and Co. 

Notation convention: The V (and below. A) symbols indicate the nabla and Laplacian 
operators with respect to x while Voj and A^; denote the nabla and Laplace-Beltrami 
operators with respect to u. Expressions of Vo; and A^ in spherical coordinates will be 
recalled later. 
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System fl2.3p -( l231) is written in a scaled form: the time and space scales have been 
chosen such that the particle speed |a;| is exactly 1 and that both d and z/ are of order 
unity. With these so-called microscopic scales, the typical time and distance between two 
particle interactions are both 0(1). We refer to Ref. [27] for a discussion of this point. 

By contrast, model (12. ip . (12. 2 p is designed to capture the large scale dynamics only, 
while averaging out the microscopic scales. Therefore, the passage from model ( 12.3p -( !231) 
to model (12.11) . (12. 2p requires a change of scales. Let e ^ 1 be a measure of the ratio 
of the microscopic length scale to the size of the observation domain. Here, the relevant 
scaling is a hydrodynamic scaling, which means that e is also equal to the ratio of the 
microscopic time scale to the macroscopic observation time. To express model (I2.3p - (l2.5p 
in terms of the macroscopic time and space scales, we perform the change of variables 
X = ex, t = et. In doing so, we must assign a scale to the interaction kernel K (i.e. to the 
interaction range R). A key assumption in the present work, following Ref. [27], is that 
this interaction range is microscopic, and therefore of order e when using the macroscopic 
scales. This change of variables leads to (dropping the tildes for clarity): 

eidtf + u ■ vr ) = -v^ ■ (F^n + rfA^r , (2.6) 

F''{x,uj,t) = u{cos 6") {Id- Lo^uj)^'', cos 6" = uj ■ u"^ , (2.7) 

x-y 



V f{y,v,t)dydv . (2. 



It is an easy matter to see that has the following expansion: 

= 9f + e^u}l + 0{e^), (2.9) 
fi^(x,t) = ^^J^, fix,t)= [ vnx,v,t)dv, (2.10) 



u}'2 = }C{id-n' ^n')^. (2.11) 

is the direction of the local flux and the constant /C depends on the interaction 
kernel K through: 

Accordingly, can be expanded: 

F' = F^ + e'F^ + 0{e'), F^ = Ff + Ff, (2.12) 
= z/(cosr)(Id-u;®w)^]^ cose' = uj-n', (2.13) 
Ff = iy{cose'){ld-uj®uj)ool Ff = u' {cos 6') {u ■ u'^) (Id - u ® uj)n' , (2.14) 

where u^cosO) is the derivative of i/(cos6') with respect to cos 6'. Because of the depen- 
dence of K upon the distance \x — y\ only, all odd powers of e vanish in the expansion. 
This would not be the case if we considered more general kernels such as those [21] de- 
pending on the angle between u and x — y . The consideration of more general kernels is 
left to future work. 
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Consequently, we will consider the following expanded Fokker-Planck model 

dtf + u-Vr + eV^ ■ {F^n = ^(-V^ ■ {FU') + rfA^r) + 0{e^), (2.15) 

where the terms Fq and F| are defined by (12.131) and (I2.14p We note that, at leading 
order, the interaction force Fq only depends on the local fiux and that the corrections 
due to the nonlocality of the interaction force appear in the 0{e) terms only. This is due 
to the assumption that the radius of the interaction region is very small (of order e) in 
the macroscopic variables. 

In Ref. [27], it has been proved that model (12. ip . (12. 2p is the formal hydrodynamic 
limit e — >^ of the mean-field model (I2.15p . Additionally, Ref. [27] provides the connection 
between the coefficients Ci, C2 and C3 of the macroscopic model to the coefficients v and d 
of the microscopic one. The goal of this paper is to investigate what diffusive corrections 
are obtained when keeping the 0{e) corrections in the Chapman- Enskog expansion of . 
These terms describe the response of the system to the appearence of gradients of the 
state variables p and VL. Here, because the fluid is anisotropic, and has only invariance 
through rotations about fi, these gradients must be split into their components parallel 
and perpendicular to Vt. 

To this aim, we denote by 

0^ = \d-VL®VL, 0\\=VL®9^, (2.16) 

the orthogonal projection matrices onto the plane normal to Vt and onto the line spanned 
by VL respectively. For a given vector X G M^, we recall that 

o^x = X - (X ■ fi)fi = ^] X (X X fi), 0\\X = (X • 

Using these projections, any vector fleld A and tensor fleld B can be decomposed into 
parallel and transverse components according to: 

A = A^ + A\\, S = 5^,^ + 5^,11 + (2.17) 

defined by 

Ai_ = O^A, A\\ = 0\\A, 

B±,± = 0±B0±, 5^,11 = 0±B0\\, = 0\\B0±, B\^ = 0\\B0\\. 

Now, we decompose gradient fields according to their parallel and normal components 
to f2. For a scalar function /, we define the normal and parallel gradients as 

V^/ = (V/)x, V||/ = (V/)||. 
Similarly, we may decompose the gradient of a vector field u into 

Vm = V_L,^M + V_L,||M + V||,_lM + V||,||U, 
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using the tensor decomposition (12.171) . Applying this decomposition to Vfi itself, we find: 



v±,±n = vn-n^{n-v)n, Vi±n = n^{n-v)n, (2. is) 

V±,||fi = 0, V||,||fi = 0. (2.19) 

The last line is a consequence of (Vf2)f2 = 0, which is found by taking the derivative of 
the relation = 1. 

In compressible Navier-Stokes equations [21], the diffusion terms can be expressed as 
functions of only two quantities constructed with the gradient of the velocity field u: the 
traceless rate of strain tensor a{u) = Vu + (Vm)^ — (2/3)(V • M)Id, and the divergence 
field V ■ u (the exponent T denotes the matrix transpose). Here, the anisotropy of the 
problem gives rise to different diffusivities in the directions parallel or normal to Q and 
we need to split the matrix Vfi into a larger number of separate entities. To this aim, we 
note that 

V ■ 1] = = Tr(V±,xfi), (2.20) 

where 'Tr' denotes the trace of a tensor. The traceless tensor V_L,_Lfi — (1/2)(V ■ Q)0± is 
decomposed in its symmetric and anti-symmetric parts a{Q) and r(r2): 

= 0±{VQ + {VQf - {V ■Q)ld)0±, (2.21) 

T{n) = = o^{vn~{vnf)o^. (2.22) 

These relations will be used in the form: 

Vx,±fi = ^{cr{n) + T{n)) + ^{v-n)o±, (2.23) 

(V^,±^^)^ = ^{a{Q)-T{n)) + ^{V-Q)0^. (2.24) 

The non-zero block of a{fl) is a 2 by 2 symmetric traceless tensor and the non-zero block 
of T{Q) is a 2 by 2 anti-symmetric tensor. We note that, for a given vector X G M^: 

T{n)x = -{T{n)fx = ((V X ^]) ■ ^]) X X n. (2.25) 

Similarly, through f l2.18p . V_Lj|fi depends only on {Q ■ V)Q and we have: 

{n-w)n = {\/nfn = xn) xn, (2.26) 

where V x ri denotes the curl of Q. Physically, {Q ■ V)Q describes the rate of tilt of Q 
as one moves along the flow lines (see flgureH]). The other quantities describe elementary 
flow patterns in the plane normal to Q: V ■ refers to convergent or divergent flows in 
the direction normal to fl while refers to swirling patterns around fl (see flgureE]) 
and a{Q) to shear patterns with one converging and one diverging orthogonal directions. 
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X 



Figure 1: (Jl • V)ri describes the rate of tilt of 17 as one moves from one point s to the 
neighbouring point s + rfs on the flow lines. In the left figure, (O • V)0 will be aligned 
with the X axis ; in the right figure, with the y axis. 




Figure 2: V ■ f2 refers to converging (or diverging) flow patterns in the plane normal 
to r2 (left picture) while r(r2) refers to swirling patterns around 17 (right picture). f2 is 
directed across the plane of the picture and pointing towards the observer (represented 
by the circled point at the origin). 




Figure 3: a{^) describes shear in the plane normal to f2. The left picture depicts the case 
associated to a\ and the right picture, the case associated to (72. is directed across the 
plane of the picture and pointing towards the observer (represented by the circled point 
at the origin). 
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Restricted to the plane normal to Q, a{Q) is a symmetric traceless 2x2 matrix. Therefore, 
it can be expressed as a linear combination of the two elementary matrices: 





M \ 




r 






K -1 


, Cr2 = I 







each corresponding to an elementary flow pattern (see figure [3]). 

In the compressible Navier-Stokes equations, there are no diffusion terms in the mass 
equation and diffusion terms in the energy equation are expressed in terms of the temper- 
ature gradient VT as a whole. Here, the temperature is constant (it is fixed by the noise 
level) but, because of the lack of momentum conservation, the diffusion terms in the mass 
conservation equation are not zero. Additionally, both the mass and velocity diffusions 
depend on Vp as well as on Vf2. Therefore, similar to Vfi, we decompose Vp into its 
parallel and normal components (respectively (fi ■ V)p and V_lp). 

Finally, the diffusion terms are composed of two parts: the first one is a quadratic 
form of the gradients in the set {V_lp, {fl ■ V)p, (fi • r(i7), V ■ Q} with coeffi- 

cients depending on (p, Q) ; the second one is a linear combination of second derivatives 
constructed by taking parallel derivatives i7 ■ V or perpendicular derivatives V_l of the 
gradients in the above list. This is precisely stated in the following theorem, which con- 
stitutes the main result of this paper: 

Theorem 2.1 (formal) The following model 

dtp + V -{cipn) =eRu (2.27) 
p {dtn + C2{n ■ + cs V±p = eR2, (2.28) 

where Ri and R2 given below, provides a second order approximation of the moments of 
the solution p^ and Q'^ of the initial model ^2.15) . The right-hand- sides are given by: 

Ri = (3v- {{n ■ Vp)n) + 7V ■ (p(v ■ (2.29) 

R2 = Q + V, (2.30) 
Q= +Qiin- Vp) Vxp + Q2 (V ■ fi) Vxp + Q3 a{n)V±p 
+ Q4 r(n)Vxp + QU^- Vp)(l^ • v)n + Qe (V ■ ■ V)Q 

+ Q7 a{n){n ■ v)n + Qs T{n){n ■ (2.31) 
V = ViO^{n-v)v^p + V2 0^{n-v){{n-v)n) + v^v^iv-n) 

+ P4 o±v ■ cx{n) + V5 o±v ■ r(n), (2.32) 

where (3, 7, Qj for j = 1, . . . , 8, Vj for j = 1, . . . , 5 are coefficients, possibly depending 
on p, which are given below. Additionally, we have (3 > 0. 

The structure of R2 is as announced: it is decomposed into a term Q which is a 
quadratic function of the gradients and a term V which consists of derivatives of these 
gradients. Both terms involve coefficients which may depend on p. The quadratic part 
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combines products of the parallel gradient of p, ■ Vp, with the perpendicular gradient 
of Pi V_lp, and similarly for ^2 (the parallel gradient of VL being VL ■ VO, the perpendicular 
ones being defined as any of gradients in the list {V ■ fi, criVt), r(f2)}) and products of the 
perpendicular gradient of p with the perpendicular gradients of VL or parallel gradients 
of p with parallel gradients of VL. The diffusive part involves only parallel gradients of 
the parallel gradient of fi, or perpendicular gradients of perpendicular gradients of f2, 
and finally a parallel gradient of a perpendicular gradient of p. In spite of its complex 
expression, R2 has a lot of structure, since a general function of this form would have 21 
different terms in the quadratic part and 10 in the diffusion part, instead of respectively 
8 and 5. 

The main objective of this paper is the proof of this theorem. In future work, the 
properties of this system will be analyzed. In particular, the question of the well-posedness 
of the system under the sole condition /3 > will be investigated, at least on a simpler 
model. Indeed, the other diffusion coefficients have no definite signs, but the constraint 

= 1 should prevent the formation of singularities, by contrast to the usual backwards 
heat equation. 

The organization of the proof is as follows. In section [3], we recall the main properties 
of the collision operator Q that are proved in Ref. [27] and provide additional prop- 
erties of the linearized operator Q about an equilibrium. This prepares the terrain for 
the Chapman- Enskog expansion which is performed in section |H We start section H] by 
an exposition of the main steps to be accomplished. Then, we successively examine the 
solvability condition for the existence of the first order correction to the equilibrium, the 
finding of an analytical expression of it as a function of elementary solutions of the lin- 
earized collision operator, and finally, the computations of the moments of this correction 
which precisely give rise to the expressions of Ri and i?2- We begin with the properties 
of Q in the next section. 



3 Properties of the collision operator and of its lin- 
earization 

3.1 Preliminaries 

We recall the expressions of the gradient and divergence operator on the sphere. Let 
X = {xi,X2,X3) be a cartesian coordinate system associated with an orthonormal basis 
(61,62,63) and let (6*, 0) be a spherical coordinate system associated with this basis, i.e. 
xi = sin ^ cos 0, X2 = sin ^ sin 0, X3 = cos6'. Let also (60, 6,/,) be the local basis associated 
with the spherical coordinate system ; the vectors ee and Ccf, have the following coordinates 
in the cartesian basis: eg = (cos 6* cos 0, cos 6* sin 0, — sin 6^), e^f, = (— sin0, cos0, 0). Let 
f{u) be a scalar function and A = AgCg + A^e^ be a tangent vector field. Then: 

Vo;/ = def ee + d^f e^, V^- A = -^dg{Ae sin 9) + -^d^A^. 
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denotes the Laplace-Belltrami operator on the sphere: 
A./ = ■ V^/ = ^deisinedef) + 



sin 9 sin 9 

We write Ffs for Fq. We introduce the 'colhsion' operator, which corresponds to the 
leading order term of fl2.15p : 

Q{f) = -V^-{Fff) + dA^f, (3.1) 
Ff = iy{ld-uj®uj)nf, (3.2) 

r2/ = "^.'^ , and j/ = / u f duj . (3.3) 

I J/ I iwG§2 

We note that Q(/) is a non hnear operator. From now on, we assume that / is as smooth 
and integrable as necessary. We note that Q acts as an operator on functions of u only and 
that the possible dependence of these functions on (x, t) can be ignored. The properties 
of Q have been demonstrated in Ref. [27] are developed in the next section. 



3.2 Properties of Q 
3.2.1 Null-space of Q 

For VL G S^, let ji = cos9 = {u ■ Q). We denote by a{^) an antiderivative of i^(/i), i.e. 
{d(7/dfi){fi) = z/(/i). We define 

,cr(tj ■ Q) 



Mn{uj)=Cexp{ 



d 



Mn(co)duj = 1. 



(3.4) 



The constant C is set by the normalization condition (second equality of (13. 4p ) ; it depends 
only on d and on the function a but not on Q. We note that, when u is constant, 
a{cos9) = V cos 9 and that Mq is the so-called Von-Mises distribution. The Von-Mises 
distribution extends the notion of Gaussian for functions defined on the sphere and is 
also known as the circular Gaussian. In the present case, the Von-Mises distribution is 
centered at (or peaked at) fi. 

The following lemma states what are the elements of the null-space of Q, i.e. what 
are the equilibria of the problem (see Ref. for the proof): 



Lemma 3.1 (i) The operator Q can be written as 

QU) = dV. 

and we have 



f 



(3.5) 



H{f) :-- 



Qif) 



Mo 



■ du = —d 



e§2 



Mo 



du < 0. 



(3.6) 
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(ii) The equilibria, i.e. the functions f{uj) such that Q{f) = form a three-dimensional 
manifold E given by 

S = {pMn{uj) I peR+, QeS^}, (3.7) 
and p is the total mass while Q is the direction of the flux of pMq{uj), i.e. 

p Mn{uj) du = p, (3.8) 



,2 



n = , JpMo = / pMn{u)udu. (3.9) 

Furthermore, H{f) = if and only if f & £ . 

An elementary computation shows that the flux can be written 

jpMn = Cipn, Ci = (cos 6')mo, (3.10) 

where for any function 5f(cos6'), the symbol {g{cos9)) m denotes the average of g over the 
probability distribution M^, i.e. 

f , . , , frT 9 (cos 6^) exp(^^^^^^) sinO dO 

{g{cose))M = / Mn{u)g{u ■ Q) duj = '[^ ■ (^-H) 

J Jq exp( ^ ^ ' ) smOdO 

fl3.10p defines the constant ci appearing in (12. ip . 



3.2.2 Generalized collision invariants 

The second set of lemmas state what are the generalized collision invariants of /. Indeed, 
we recall that the collision invariants are classically defined as the functions ip{uj) such 
that 

/ Q{f)ijdco = 0, V/. (3.12) 

However, it is readily seen that the linear vector space of collision invariants is of dimension 
one, while the hydrodynamic limit requires that its dimension be equal to the dimension 
of S, which is 3 in the present case. To find the missing collision invariants, we slightly 
weaken the definition. We fix i7 G arbitrarily, and we define a Generalized Collision 
Invariant (or GCI) associated to as a function ip which satisfies fl3.12p only for functions 
/ with direction Qf = Q. This constraint is linear and can be resolved by the introduction 
of a Lagrange multiplier. This leads to the following definition: 

Definition 3.2 Let f2 e 6e given, ipioj) is a Generalized Collision Invariant (or GCI) 
associated to Q if and only if 

/ Q{f)tljdu = 0, V/ such that Qf = Q, (3.13) 
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It is shown in Ref. [27] that, using (13 .Sp and Green's formula, f l3.13p leads to the 
following problem defining the GCI's associated to a direction fl: 3(3 G M^, such that 
/3 ■ = and 

V, ■ (MnV^V') = /3 ■ (^^ X tu)Mn. (3.14) 

This problem is obviously linear, so that the set Cn of GCI's associated to is a linear 
vector space. In a cartesian basis (ei, 62, fl) and the associated spherical coordinates {9, 0), 
we have /3-(f2 xw) = (— /3i sin0+/32 cos 0) sin 6* with = (3-ek, k = 1,2. Therefore, we can 
successively solve for ipi and ip2, the solutions of (I3.14p with right-hand sides respectively 
equal to — sin^sin^M^ and cos^sin^M^. The following lemma provides the framework 
for solving (I3.14p . It is based on Lax-Milgram theorem and is proved in Ref. [27] : 

Lemma 3.3 Let x £ i^^(S^) such that J xd^ = 0- The problem 

■ (Mj,V^^) = X, (3.15) 

o 

has a unique weak solution in the space if^(§^), the quotient of the space H^(S'^) by the 
space spanned by the constant functions, endowed with the quotient norm. 

So, to each of the right-hand sides xi = — sin^sin^M^ or X2 = cos sin ^M^ which 
have zero average on the sphere, there exist solutions ipi and ip2 respectively (unique up 
to constants) of problem (13.150 . We single out unique solutions by requesting that ipi 
and -02 have zero average on the sphere: J xpkdu = 0, k = 1,2. Then, we have, as a 
consequence of Lemma 13.31 

Proposition 3.4 The set Cq of generalized collisional invariants associated with the vec- 
tor Q which belong to H^iS"^) is a three dimensional vector space Cq = Span{l,ipi,4'2} ■ 

More explicit forms for ipi and ^2 can be found. By expanding in Fourier series with 
respect to 0, we easily see that 

■01 = —g{cos9) sin0, 02 = g{cos9) cos0, (3.16) 

where g{fi) is a solution of the elliptic problem on [—1, 1]: 

-(1 - /i2)9^(e"('^)/''(l - iJ^)d^g) + e"('^)/'^(7 = -(1 - /i2)3/2e-(^)/'^. (3.17) 

To solve this problem, we apply the following lemma: 

Lemma 3.5 Let X = | (1 - j^^y^/^g G L\-l, 1)}, V = {g e X\{1- ^^f^d^g G 
L^(— 1, 1)}. Let a{fi) belong to L°°(— 1, 1) such that there exists ao > and a(yu) > cto- 
Then, for any / G X, there exists a unique solution g eV of the problem 

-(1 - /i2)9^(e'^(^')/'^(l - i,'')d,g) + a{ii)g = f. (3.18) 

Additionally, the maximum principle holds: if f is non-positive (respectively non-negative), 
then, so is g. 
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The next lemma, in the spirit of the previous one, will prove useful in the sequel: 
Lemma 3.6 Let X = L^-l, 1), V = {g e fi^y^^df^g G L\-l, 1)}, k= {/ G 

-j^ o o o 

A* I / d/i = 0}, V= V/M. Then, for any f EX, there exists a unique solution g GV 
of the problem 

-d^{e'^^^y\l-li'')d,g) = f. (3.19) 



Both lemma are direct consequences of Lax-Milgram's theorem. The first one is proved 
in Ref. [27]. For the second Lemma, thanks to a Poincare inequality we note that the 

o 

semi-norm of V is equivalent to the norm of V on the quotient V- For both problems, no 
boundary conditions at ±1 need to be prescribed. This is due to the degeneracy of the 
elliptic operator at these points. 

Applying Lemma [331 shows that the function g, solution of fl3.17p . is uniquely defined 
in the space V. For convenience, we introduce = (1 — jJ^)~'^^'^g G 1,1) or 

equivalently h{cos9) = g{cos9)/ sin 9. We then define 

iIj^{uj) = {^l X uj) h{fi) = ipiCi + ilj2e2, jj = {uj ■ ^l) . (3.20) 

ipQ is the vector generalized coUisional invariant associated to the direction Q and is 
uniquely defined by the problem 

■ (Mf, V^4) = (^ X tu)Mn, I i^ndco = 0. (3.21) 

We note that, by the maximum principle, h < 0. 



3.3 Linearization about the equilibrium state 

We introduce a macro- micro decomposition of /: 

f = PfMn,+v, (3.22) 
Pf and Qf being the density and mean velocity direction of /, i.e. 

Pf= [ fdoj, (3.23) 

and VLf given by (13. 3p . These definitions of p/ and VLj are equivalent to saying that ip 
belongs to the space: 

= |(^ G L^(§^) I j Lpdco = and rif X j (^a;rfcj = o|. (3.24) 
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For any given G we define the following linear operator operating on ^q: 



Ma 



(3.25) 



Inserting fl3:22|l into (ES]) yields that 

Qif) = Ln^icp). 

We now precise the functional setting. Let Xq = {lp\ f {ipl"^ M^^ du < oo}, = 

{^p eXn \ J Mn\VUM^^^p)\^du < oo}, Xn= Xn/(MnM), Yn= Vn/(MnM). We 
perform the usual identification of Hq with its dual X'^. The following lemma states the 
properties of L^. 

Lemma 3.7 Consider Lq as an operator from Yn into its dualY'Q, defined by the bilinear 
form on Y^: 

{L.^, X),vi,..„, = -d M„V„ (^) ■ V„ d.. (3.26) 

Then, 

(i) The null-space of Lq, is the linear space spanned by Mq. 

(a) Let ( G Xn- The equation LQip = ( has a solution ip in Yq if and only if ( satisfies 
the solvability condition: 

Cduj = 0. (3.27) 



And Lp is unique if it additionally satisfies ^3.21^ . This unique solution is written = Lq^(^ 

and L'^ is called the pseudo-inverse of Lq_. 

(Hi) The solution (p belongs to if and only if ( satisfies 

C^nduj = 0, (3.28) 
where if^ is the GCI Ii3.20\) . 

Proof: Statement (i) is obvious. To prove (ii), we note that the bilinear form fl3.26p 

o 

defines a bilinear form on the quotient space Yn- The solvability condition fl3.27p is the 

o 

necessary and sufficient condition for an element of X^ to belong to the dual space (Xn)'. 
Then, the proof of (ii) follows from the application of Lax-Milgram theorem. Indeed, by 

Poincare's lemma, the bilinear form f l3.26p defines a norm on Yn which is equivalent to 

o 

the norm of V^. Therefore, there exists a unique solution (p GVn of the problem 
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which uniquely defines a solution in Vo provided the cancellation condition fl3.27p is 
imposed on ^p. Finally, inserting x = MqiPq in f l3.29p and using (13.211) . we deduce that 

,2 7§2 \Mn/ 

= d ■ (MnVujifn) du = d Q x (fudu. 

JS2 Mq J^2 

Property (iii) is an immediate consequence of this identity. ■ 



3.4 Coefficients of the hydrodynamic model 

We finish these preliminaries by recalling the expressions of the coefficients Ci, C2 and C3 
as they were derived in Ref. ^7]. For any functions g{cos6), h{cos6) with h > 0, we 
denote by {g)h the average of g over the probability distribution defined by h: 

J g{uj ■ n) h{uj ■ n) du f^^ g {cos 9) h{cos 6) sin 9 d9 

{g{cos9))h = jr— -— . . 3.30 

j h{uj-\l)duj Jg h[cos9) sm 9 d9 

Then, the constants ci, C2 and C3 of model (12. ip . (12. 2 p are given by: 

Cl = (cos6')Af , C2 = (cos6')(sin2 0)^M,/ , C3 = d(-) . (3.31) 



(sin2 e)uhM 



With simple computations, one can check that these definitions are equivalently expressed 
by the following relations: 



{{u-n)~ci)Mndu = 0, (3.32) 

[ ^{{uj-Q)-C2){l-{co-nf)hMQdu = 0, (3.33) 
Js2 d 

/ ~ " ■ hMnduj = 0. (3.34) 

4 The Chapman-Enskog expansion 
4.1 Setting up the expansion 

We introduce the macro-micro decomposition (13.221) in a scaled form: 

r=p'Mns+eG', (4.1) 

where p'^ = pfe and Vt^ = Vtfe are the density and velocity direction of the solution of 
the kinetic model (I2.15P and where G This leads to: 

{dt + u- V){p'Mn^ + eG') + eV^ ■ (F^p'Mn^) = L^XG') + 0{e^). (4.2) 
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This justifies the scahng (14.11) because G"^ appears as an 0(1) quantity, showing that the 
correction to equihbrium eC^ is 0{e). 

Integrating (14. 2 p and using Lemma [3.71 (ii). we find: 



§2 



(dt + u- V)(p"Mf,. + eC) du = 0{e'). (4.3) 
The contribution of the term p^M^e has been computed in Ref. [27] and we get: 

dtp'' + V ■{cip''Vt')=eR{ + 0{e^), R\ = -V-l G'todco. (4.4) 
Indeed, since G"^ G we have 

dtG' doj = dt{ I G'duj] =0. 



Now, multiplying (14.21) by the GCI ipn{^) = (f^ x w) h{Q ■ u), integrating with respect 
to u and using Lemma 13.71 (iii) and that G , we find: 



X [{dt + CO ■ V){p'Mn^ + sG') + eV^ ■ {F^p'M^.)] hudu^ = 0{£' 



^ (4.5) 
Now, thanks to the computations of Ref. [2^, we have: 

fi" X I (9i + W ■ V)p"Mf,. hujdu^^= ^ ^Q^^e ^ 

+C2(r^' ■ v)fi") + C3( Id - fi= ® ^]")Vp"} , 

from which we deduce that 

(9ffi" + C2(l^" ■ V)fi") + C3( Id n')Vp' = eR2 + 0(5^), (4.6) 

2d r 

R2 = - , . o , , , / [dtG' + uj-VG'' + 

(sm 9uh}Mae J§2 

+V,-(F2VMnO]/ic^rfc^. (4.7) 

If we omit the 0(5^) remainders in eqs (14.41) and (14.71) . we find a macroscopic model 
which approximates the moments of the Fokker-Planck model (I2.15P up to 0(5^) terms. 
The goal is now to compute Ri and R2 and to show that, up to 0{e) terms, they have 
the expressions given by Theorem 12.11 Obviously, this requires the computation of G^. 

From (14. 2p . we find that = + 0{e) where G^ is a solution of the problem 

LnAGn = {dt + uj-V){p'Mns). (4.8) 

Therefore, replacing G'^ by G*^ in the expressions of Ri and R2 will not change the order of 
the approximation. From now on, we will omit the tildes and consider G^ as the solution 
of g3). 

Now, the plan is as follows 
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1. Show that (14.81) has a unique solution belonging to the space This means 
proving that the right-hand side of (14. 8 p satisfies the solvability conditions (I3.27P 
and dSSHl). 

2. Compute the expression of G^, i.e. invert (14. 8p . 

3. Insert the expression of in the definitions of -Ri and R2 and compute them. 

We now successively perform these tasks. In the following, we will omit the exponents e 
to make the expressions lighter. 

4.2 Preliminary lemmas 

We will have to deal with the integrals over uj of expressions involving tensor products of 
Cj_ and 0\\ . To this aim, we use a spherical coordinate system associated with a cartesian 
basis whose third basis vector coincides with VL. We denote by {9, (p) the associated angular 
coordinates as defined in section ISTTl For a function u{u), we define 

1 f'^'^ 
{u)<f> = ^ J u{e,(f))d(f), 

the average of u over the angle 0. The proof of the following lemma is easy and omitted: 
Lemma 4.1 (i) For all odd tensor powers p, we have 

iio^Lorn<p = o- (4.9) 

(a) The first even tensor powers ofO±u are given by: 

{{0^u;f')^ = ^sm'eO^, (4.10) 
{{0^uf')^ = lsm'eO^, (4.11) 

o 

where we define the fourth order tensor 0_l .■ 

= {0,r + ii0.r),,,um + ((^-)"')[28M46] ' (4.12) 

and the subscript [ij] denotes contraction with respect to indices i and j . Restricted to the 
plane normal to Q, 0± can be written: 

(0_L)ijfcz = dijdki + SikSji + 5ii5jk, (4-13) 

where ((5ij)i,je{i,2} is the two-dimensional Kronecker tensor. We note that 0± is invariant 
under all rotations of the plane normal to Q, i.e. it satisfies 

{0±)ijkiRii'Rjj'Rkk'Ru' = {^i.)i'j'k'i', 

for all rotations Rw of the plane normal to Q, where Einstein's repeated index summation 
rule is assumed. 
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4.3 Sovability of the equation for G 

We state the first lemma: 



Lemma 4.2 We have: 



{dt + uj-V){pMn) = A|| ■ V||p + A^- VxP 

+ p(5x,±: Vx,±^^ + 5||,^: V||,x^^) + 0(£), (4.14) 



where 



A^ = Mn{l-'^)o^uj, (4.15) 

= Mniiuj ■ n) - ci)n , (4.16) 

= Mn C^iOLUj) ® {O^u) - CiOx) , (4.17) 

Bi± = Mn^iioo ■ n) - C2)n ® {0±uj) , (4.18) 
and where ':' denotes the contracted product of two tensors. 
Proof: We have 

{dt + u- V)(pMo) = Mn (^{dt + u-V)p + P ^^^l^^'h dt + oj ■ V)n 

Classically, in the Chapman-Enskog procedure, time derivatives are replaced by space 
derivatives, using the following identities 

dtp = -V-{cipn) + 0{e), 
pdtVl = -C2p{n ■ V)n - C3 OxVp + 0{e), 

which are deduced from (14.41) and (14. 6p . For any tangent vector Q to §^ at Q, we have: 
Then, we note that: 

a; ■ Vp = 0±UJ ■ V_lP + {OO ■ fi)^^ ■ V||P, 

v-n = o^: (Vx,±fi), 

uj-{n- V)fi = (fi O 0±uj) : V||,xfi, 

u-{uj- V)Q = {O^u ® O^u) : Vx,xfi + {uj-Q){Q® O^u) : V||,x^- (4.19) 

Collecting these identities, we find expressions f l4.15p to (14.181) . ■ 
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Lemma 4.3 The quantities A^, A\\ B±^±, B\\^± satisfy (separately) conditions \3.21\l and 



^3.2^) . As vectors or tensors, this means that they satisfy these conditions componentwise. 



Proof: We summarize the main arguments and leave the computational details to the 
reader. 

(i) A_[_ satisfies (I3.27P because of ( 14. 9p and (13.28^ as a consequence of (]3.34p . 

(ii) A\\ satisfies f l3.27p as a consequence of fl3.32p and f l3.28p because of (14. 9p . 

(iii) -B±,± satisfies (I3.27P as a consequence of (I3.32p (after an integration by parts with 
respect to 6) and (I3.28P because of (14.90 . 

(iv) X satisfies (I3.27P because of (14.90 and (I3.28P as a consequence of (13.330 . ■ 



4.4 Computation of G 

We have shown that the right-hand side of (14.80 can be decomposed into four different 
terms, corresponding to derivatives of p and Vl in the directions normal or parallel to f2, 
and that each of these four components satisfies the solvability conditions (I3.27P and (I3.28P 
separately. We now compute the pseudo inverse L'^ applied to these four components. 

Lemma 4.4 We have: 

Ai. := -L^^A^ = Mn a± 0±u , (4.20) 

A\i := -L^^A\\= Mna\\n, (4.21) 

B±,± := -i^n'5±,± = Mn {h (O^co) ® (O^co) + h O^} , (4.22) 

Bi± := = Mn fey n ® (O^u) , (4.23) 

where a±, a\\, hi, 62 ond h\\ are functions of u ■ Q. They are defined by the following 
relations (letting p = u ■ Q): 

1. a± = a±{p)^yi — p? is the unique solution of ^3.13) with 

a{p) = et, /(/.) = [l - ^) (1 - p-fl\ (4.24) 

2. a\\ is the unique (up to an additive constant) solution of 1^3. 19i) with 

/(/i) = ie5(/i-ci), (4.25) 

and the constant is adjusted in such a way that J^2 A\\ du = 0. 

3. 61 = 61(1 — /i^) is the unique solution of Ii3.18\) with 

a(/i)=4ef, /(/.) = ^e^(l-/i^)2. (4.26) 
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4- 1)2 is the unique (up to an additive constant) solution of ^3. 1 9) with 

/(/.) = et(26i-|), (4.27) 

and the constant is adjusted in such a way that J^^ B± du = 0. 
— /i^ is the unique solution of Ii3.18\) with 

«(/^) = e^, /(/^) = ^e5(^-C2)(l-/iT/'- (4.28) 



Proof: Preliminaries: using spherical coordinates, we check that if ip is of the form 

(f = MnCfc(cos6') cos k(f), then, 

-Lnip = ——^ cos kct) [-(1 - ii^)d,{e''l\l - fi')d,Ck) + ^e'^l'^Cu] , (4.29) 



with n = cos 9. Similarly, ii ip = Sk{cos9) sin kcj), then Sk satisfies the same identity 
with cos k(f) replaced by sin k(f). 

Proof of (i): It is a matter of computation to show that A_i_ defined by (14.201) is 
a solution of —LfiA± = A± provided that d± satisfies (13.181) with data a and / given 
by (14.241) (use (I4.29P with k = 1). Now, it is clear that A± satisfies the normalization 
condition (13.271) . because of (14.91) . Therefore, A± is the unique solution called —L^^A±. 

(ii) Using (14.291) with k = 0,we show that A\\ defined by (I4.2ip is a solution of -LqA\\ = A\\ 
provided that ay satisfies (I3.19P with data / given by (14.250 . ay is defined up to an additive 
constant, which means that A\\ is defined up to the addition of a function proportional to 
Mu^. The coefficient can be chosen in such a way that condition (I3.27P is satisfied. This 
solution is the unique solution called — L^^Ay. 

(iii) We proceed similarly. Using (14.291) successively with k = 2 and = we show 
that —LqBj_^± = -B±,_L provided that bi and 62 are specified as stated in the theorem. 
Additionally, with (I4.10p . we have 

(^x,x)0 = (6i^ + 62)O±, 

and since 62 is defined up to a constant, we can adjust this constant to satisfy the nor- 
malization condition (I3.27p . The so-defined is the unique —L^^B±^^. Note that, 
because of the factor 1 — /x^ in the expression of / in (I4.26p . it is an easy matter to show 
that bi = 61/(1 — /i^) belongs to and that the assumptions for the application of lemma 
13.61 are satisfied. 

(iv) We proceed exactly in the same way for Using (I4.29P with A; = 1, we find that 
B\\^± is a solution of —LqB\\^± = provided that b\\ satisfies (I4.28p . The normalization 
condition (I3.27P is satisfied because of (14. 9 p which proves that the so-defined B\\^^ is the 
unique —Lq^B\\^±. u 
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Lemma 4.5 The following relations are satisfied: 



{a^ sin^ 9)Mn = 0, (a||)M„ = 0, (4.30) 
{^h sin^ 9 + b2)Mn = 0, ik sin2 9)Mn = 0. (4.31) 



Proof: Since A±, A\\, . . .belong to the space their integral against 1 and x u over 
a; e vanishes (see definition (13.241) ). Using lemma WA] this leads to the above listed 
relations. ■ 

Finally, as a consequence of lemma 14.41 we can summarize: 

Lemma 4.6 we have: 

-G = -L-^\{dt + uj-V){pM^)) 
= i|| ■ V||p + i± ■ Vxp + p(5x,x : V±,±l] + 5|,,x : V||,x^^) +0(£). (4.32) 

We can decompose G into even and odd powers of 0±u and write 

G = Ge + Go + 0{e), (4.33) 

-Ge = i|| ■ V||p + p5x,± : v±,i_n 

= Mn {a|, n-Vp + p[bi {{O^uj) ® (O^uj)) : Vfi + 62 : Vfi] } , (4.34) 
-Go = i± ■ V±p + pB\\^i_ : Vii_Q 

= Mn {a± 0±uj ■ Vp + pfoy (n ® (O^tu)) : Vfi} . (4.35) 



In the sequel, we will omit to mention the 0{e) remainder. It should be understood that 
all results are up to a term of this order. 

4.5 Computation of the 0(e) corrections 
4.5.1 Computation of Ri (H^D 

In this section, we compute Ri, the right-hand side of (14. 4p . Its expression is given in the 
following statement: 

Lemma 4.7 Ri is given by formula \2. 29\] with 

X 

13 = {a\\ cos 9) 1 = {{-bisiv? 9 + b2)cos9)Mn- (4-36) 
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Note: compared with fl4.30p and f l4.3ip . there is an additional factor cos^ inside the 
brackets. 

Proof: We multiply f l4.32p by uj and integrate over uj. But, because G G $n; the normal 
component of J^2Guduj to Q vanishes and the projection upon Q is the only non-zero 
component. It is obtained by multiplying f l4.32p by {Q ■ and integrating upon u. In 
this integration, the contribution of the odd part Go vanishes by (14.90 . The contribution 
of the even part Ge is readily found to be f3{Q ■ V\\p)Q + 'yp{0± : (Vfi)_L,_L) Q. The proof 
is ended by using f l4.19p . ■ 



Lemma 4.8 We have /3 > 0. 

Proof: We can write: 

P = I a\\ Mn {oj ■ Q) du = / a\\ Mn {{u ■ Q) - ci) du 



ys2 



y§2 



:a|| Mn) Ln(a|| Mn) M^^ du > 0, 



by the non-positivity of (see (I3.26P ). In the second equality, we have used that a\\ Mq 
satisfies (I3.27P (see also (I4.30p ). The third equality is obvious and the fourth one is just 
using the definition of ay Mq (see fl4.2ip ). (3 is strictly positive, otherwise, a\\ would 
belong to the kernel which of Lf^, which is spanned by Mq. Since, besides that, a\\ Mq 
satisfies (I3.27p . it would be identically zero. But, applying Lq to it, then {{u ■ fl) — Ci) Mq 
would also be identically zero, which is obviously not the case. This concludes the proof 
by contradition. ■ 



4.5.2 Computation of R2 diTfl) 

Lemma 4.9 R2 is given by Ii2.30\) . The coefficients are given by (J^J^j ( |^.50[ j 



Proof: We first compute the term involving dfG: 

T = -O^ [ dtGhudu. 

J§2 

We have 



T = -0±dt(^j^ Ghudu^+0± (^J^ Gh'u^udu^ dt^ = T1 + T2, 
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where h' denotes the derivative of h with respect to ^ = u ■ VL. For Ti, we decompose u 
according to transverse and normal components: 

Ti = -O^ dt (^j^ Go h O^co du?j - 0± dt (^j^ Geh{uj-n)n du?j = + T^, 

where we have used fl4.33p and f l4.9p to introduce the even and odd parts of G. Thanks 
to ( KTO\i and flCTj) . we find 

with 

= sin^ ^ ^)Mn, = sin^ ^ ^|| h)Mn, 
Xli = (cos 6* ay h)Mn^ ^12 = sin^ 6* COS 6* 61 + cos 6* 62) h)MQ- 
We proceed similarly for T2. Since u ■ dtfl = {0±u) ■ dt^, we find: 

T2 = (^J^ Ge h' {O^u ® O^uj) duj^ dtQ. 

And thanks to (00]) . KIT} and f Hl3|) . we find 

T2 = -A21 {n ■ Vp)dtn - A22P(v ■ n)dtn - X2spain)dtn, (4.37) 

with 

1 11 

A21 = (-sin^6'a|| h')Mn, A22 = ((-sin'' 6* 61 + - sin^ 6* 62) /i')Afn, 

A23 = (^ sin'' 6161 h')Mn, 

o 

where we have used that 

= {a{n) + 2(V ■ n)0±),j . (4.38) 
Now, we note the following relations: 

o±idtO^\/p) = o^dtVp - (n ■ vp)dtn , 

o±{dt{p{n ■ v)fi)) = dtp{n ■ v)fi + p{dtn ■ v)n + pO±{n ■ v)dtn , 

= ^{a{Q) - r{Q))dtn + ^(v ■ Q)dtn , 

o±dt{{n ■ Vp)fi) = {n ■ Vp)dtn , 
o±dt{p{v ■ = p(v ■ n)dtn . 
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Collecting all these identities, we get: 

T = a; o± vdtp + A2 {n ■ vp) Otn + A3 {{n ■ dtp + \'^pa{p)dtn 

+ \'^pV{n)dtVt + x'^piv ■n)dtVL + \'TpO^{yi-v)dtVL (4.39) 

with 

Now, we use that up to order e terms, we have: 

dtp = -V ■ (cipfi), 

pdtQ = -C2p {fl ■ V)n - C3 (Id - (g) fi)Vp, 
and replace the time derivatives appearing in fl4.39l) by space derivatives. Using that 

o^{vn)Vp = (Vfi)x,±Vxp = ^Hn) + T{n))v±p + ^(v ■ n)v±p , (4.4o) 

o^{n ■ v)Vp = o±{n ■ v)v±p + {n ■ Vp){n ■ v)n , (4.4i) 

pO^in ■ V)^ = ■ V) V±p - ■ Vp)Vxp . (4.42) 
We then get: 

T = x'l (V ■ n)v^p + A'2'pVx(v ■ n) + x'^a{n)v^p 
+ x'iT{n)v±p + x'^o±{n-v)v±p + x'^{n-Vp){n-v)n 

+ Xj-{n-Vp)v±p + x'^piv ■n){n-v)n + x'^pa{Q){n-v)n 

+ A;'oP^(fi)(^] ■ v)n + x';^pO±{n ■ ■ v)n) , (4.43) 



with: 



■^1 — ~'^i~2 "^6*^^ ' ^2 — "-^'1^1 , A3 — "-^1"^ — -^4^3 5 

A4 = "-^1"^ -^5^3 ) -^5 = —^[Ci — A7C3 , Ag = —X[ci — A2C2 — A3C1 

K = -^2^3 + A7C3 , Ag = -A3C1 - AgC2 , Ag = -A4C2 , 

Now, we turn towards the term involving u ■ VG: 

S = -0± [ {u-V)Ghudu. 

JS2 
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We decompose 



5" = -O^V-(^j G hco ^uduj + 0± (^j G h' u ^ co ^ co dco^ (VQ) 



— Si + S2, 

where the second term has the following meaning: 

G h' tu (^u ®u 



duj] (Vn)] =( / Gh'uj®uj®ujduj] (yn)jk, 
J ) i Vis2 ) .-^-^ 

and Einstein's summation convention is assumed. We again decompose 5*1 into parallel 
and normal components: 

Si = -Oi_V ■ (^j^ Geh (O^co) ® (O^co) du?j - 0±V ■ (^j^ Gah {u ■ Q) Q ^ {0±uj) du^ 
-C± V ■ (^J^ Goh{uj- Q) {0±uj) ^Qdu?j - 0±V ■ (^j^ G^h {u ■ Qf Q ^ Q du?j 

= sl + ... + st, 

where again, we have used (14. 9 p to restrict to the even (Ge) or odd (Go) components of 
G with respect to (O^u). Using similar computations as for f l4.38p and ( I4.37p . we find 

- f Geh {O^u) ® {O^u) du = r/i\ {Q ■ Vp) + 7]\^ p (V ■ fi) + P^(fi), 

JS2 



with 

= (l sin^ e a,, h)M,, VI2 = {\ sin^ Obih + l sin^ 6 h /.)m„, 
vis = (gSin^^fei /l)Mn- 

We now note that 

CxV ■ {{n ■ Vp)Ox) = c±(vn)Vp + o^{n ■ v)Vp - {n ■ Vp){n ■ v)fi, 
c±v ■ (p(v ■ n)o±) = v±p(v ■ n) + pv±(v ■ n) - p(v v)fi, 

CxV ■ (p(T(fi)) = a{n)V±p + pCx(V • a{n)). (AAA) 
But then, with fICT]) and (gZH), flCTj) gives 

Collecting the above identities, we find 

^i' = (^^n + ^L) o^(^^)Vxp + \]\i r(fi) Vxp + i^-n\i + (v ■ ^]) v^p 

+r/i\ ■ V)VxP + r/ls pVx(V ■ ^) - r]l^ p(V • • V)n + r/jg pOx(V ■ 



25 



Now, turning to Sf we have: 

Goh{uj ■ {0±uj) doj = r]l^n® V±p + 77^2 P ^ ® ■ V)fi, 



with 



Noting that 



ril^ = {-siii^ 6 cos 6 a±h)Mn J ''?i2 = (2 ^^^^ ^ '"'-'^ ^ ^11 (4.45) 



c^v ■ (n ® Vxp) = (V ■ n)v±p + ■ v) v^p, 

c^v • (^] ® (f] ■ v)f]) = (V • ■ v)n + o±{n ■ ■ v)fi), 



we get 



Si = vu (V ■ n)V^p + Vn 0±{n ■ V)VxP 
For S'f , we have 

- / Goh{uj-Vt) {O^u) ®nduj = rjl^ V±p ®Vt + rjl^p {{Vt ■ V)Vt) ® Vt, 
with rjli and 77^2 given by (I4.45p . With 



= ^aim^p - ^r(fi)Vxp + ^(v • m±p, 



we find 



Then, for Sf, we write: 

- / Goh{u- doj = T]^^ {Q ■ Vp) + 77^2 P (V ■ fi), 
with rjf-^ and 77^2 given by 



T]^^ = (cos^ 9 a|| /i)Afn, ^7^2 = ((o sin^ ^ ^1 + ^2) cos^ h)Mn- 
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With 

o±v ■ (p (V ■ n)n ®n) = p(v v)fi, 

we find 

= r]^^ {n ■ Vp) + r^l^ p( v ■n){yi-v)vt. 

Now, we turn to 5*2. Using that (Vfi)r2 = 0, the decomposition of a; into 0±_uj and {(jj-VL)VL 
reduces to: 

82 = ( I Goh' 0±uj ® 0±uj (g) 0±uj du) {VQ)jk + 

\Js'^ J ijk 

+ ( [ Geh'{Lo-n) o±uj ®n® o^Lo dco] ivn)jk = si + si 

(where again, Einstein's summation convention has been used). Using (I4.38p . we find 
-Si = r]l^ {a{n)V±p + 2(V ■ fi)V±p) + r72^p(a(fi)(fi ■ V)^] + 2(V ■ • V)fi), 

-Si = vli ■ Vp)(fi ■ w)vt + ri^pa{n){{yi ■ + r]l^p{v ■n){yL- v)fi, 

with 

Till = sin^ 6 a± h')Mn, VI2 = sin^ ^ h^')M^, 

V21 = sin^ e cos 6* ay h')Mn, V22 = sin^ 6 cos 61 61 h')Mn, 

Z o 
1 1 

V23 = ((^sin^6'6i + -sin^6'62) cos9h')Mn- 
Collecting all these identities, we find 

S = r][{V- n)V±p + r]'^ pV±{V ■ + % cT{n)V±p 
+ r]'^T{n)V±p + r]'^0^{n-V)V^p + r/;(fi- Vp)(fi- V)fi 

+ r/^p(v-fi)(fi- v)fi + 7]'^pa{n){yt-v)vt + r][^pO^{yL-v){{yt-V)n) 

+ r^[^pO^V ■a{n), (4.46) 



with 



13 11 

V'l = + + ~ '^'^21 , ^2 = ^1^2 , ^3 = 2^11 + ^i^3 + 2^11 ~ ^21 , 

/-'-I -'-2 / 1,2 / 4 2 
= 2^11 - 2^1^ ' " ' " ~ ' 

^8 = + VI2 + ^^2 - 2r?22 - ^23 ) ^9 = "^22 " ^22 . V'u = "^Vu > ^12 = ^13 • 
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We now turn to the last term 



U = -O^ [ ■ {F2pMa)hLodoj := Ojj. 

JS2 



The fc-th component U ■ Ck oi U m a, Cartesian basis (6^)^=1^2,3 can be transformed by 
using Stokes theorem on the sphere: 



U-ek = - V^-{F2pMn)h{u-ek)du= / pMn ■ Vu^)ih{u ■ Ck)) du. 

An easy computation gives 

{F2 ■ V^)ih{Lo ■ Cfc)) = {cu ■ ek)h\F2 ■ n) + h{F2 ■ e^). 

Therefore, 

U = Ui + U2, 

Ui=pO^[ Mn{F2-n)h'LodLo, U2 = pO±[ MnhF2dco. (4.47) 

From fl2.12p . we can write F2 = F2 + F2. Introducing this decomposition into the expres- 
sions fl4.47p of Ui and U2, we get 



Ui = ul + ul, U2 = U] + f/. 



2 

2 1 



where for instance, is defined by the first expression (14.471) with F2 substituted by F^. 
In each of the expressions defining U'^ with p,q = 1, 2, we decompose cu into 0±oj + {uj-Q)Q 
and, using (14. 9p . keep only the even powers of 0±u. We find: 

with 

= -(^sin^^ cos6uh')Mn, = sin^ ^ ^')Mn, 
C2 = ((1 - ^sin^^)z//i)Af^, ^2 = -(^sin^^ cos9u'h)Mn, 
Therefore, using (12. lip and the fact that j = cipQ (see Ref. [27]), we have: 

Next, we decompose: 

o±A{pn) = o±{2{vnfVp + pv ■ (vn)). 
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Now, with ( 12:231) and (127781) . we have: 

o^y ■ (v^]) = ^OxV ■ aisi) + ^o^v ■ r(fi) + ^o^v ■ ((v ■ ^)o^) 

+CxV■(fi®(^]■v)^]), 
OaN ■ ((V ■ = v±(v ■ 1]) - (V • ^])(^] ■ 
o^y v)vl) = (V • ■ v)n + o^{n ■ ■ v)fi. 

Collecting these results, we finally get: 

U = ^,{V- Q)V^p + 6pVx(V ■ Q) + ^3cr{Q)V^p 

+ eiip(^±(fi- V)((fi-V)fi), +62pO±V-a(n) + 63P0±V-r(fi), (4.48) 



with 



6 — ^ ) ^2 — 2^ ' ^3 — ^ , ^4 — — ^ , ^6 — 2^ , ^8 — 2^ ' 

_ 1 _ 1 

^11=^5 ^12 = -i , 63 = 2'^ ■ 

We can now collect (l443l) . (I446l) . (j448l) and insert them into (l230l) and find 
R2 = Ci (V ■ fi)VxP + C2pVx(V ■ 9) + C3fT(fi)V^p 

+ C4r(fi)Vxp + C5C^(^^- v)Vxp + C6(^^- Vp)(fi- v)r] 

+ C7-(^^- Vp)v±p + C8P(v-i])(fi- v)fi + C9PfT(fi)(^^ ■ v)^] 
p 

+ Ciopr(^^)(fi- v)fi + CiipC±(fi- v)((r]- v)fi) + Ci2pC±v-a(fi) 
+ Ci3P0±v-r(fi), 

with 

and where we have defined the missing coefficients A" for j = 12, 13, rjj for j = 7, 10, 13 
and for j = 5, 7, 9, 10 as zero. 

Now, the proof is complete. The expressions of the coefficients Qi and Vi are as 
follows: 

Qi = -Ct, Q2 = Ci, Qs = C3, Q4 = C4, 
p 

Q5 = Ce, Qe = PCs, Qr = pCq, Qs = pCio, (4.49) 
^^1 = C5, 1^2 = pCii, '^3 = PC2, = pCi2, = pCi3- (4.50) 
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5 Conclusion 



In this paper, we have provided the 0{e) corrections to the hydrodynamic model derived 
in Refs. [26| |27] from the kinetic description of the Vicsek ahgnment dynamics [53]. The 
0{e) corrected model involves diffusion terms in both the mass and velocity equations as 
well as terms which are quadratic functions of the first order derivatives of the density and 
velocity. To express these terms, it is necessary to decompose the density p and velocity Q 
and their gradients in the directions parallel and normal to Q, thereby expressing that the 
fluid is non-isotropic about Q. Future works are concerned with the mathematical theory 
of this system at least in a simplified form, with the derivation of asymptotic formula 
for the coefficients in the limits of small and large noise and with numerical simulations 
and comparisons with the original particle dynamics. In particular, a question to be 
examined is whether including the 0{e) corrections in the simulation allows to bypass the 
ambiguities of the non-conservative hydrodynamic model (see Ref. and to yield a 

better approximation of the solutions of the original particle model 
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